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ABSTRACT 

Ambipolar diffusion (AD) is believed to be a crucial process for redistributing magnetic flux in the 
dense molecular gas that occurs in regions of star formation. We carry out numerical simulations of 
this process in regions of low ionization using the heavy ion approximation. The simulations are for 
regions of strong field (plasma /3 = 0.1) and mildly supersonic turbulence (A^ = 3, corresponding 
to an Alfven mach number of 0.67). The velocity power spectrum of the neutral gas changes from 
an Iroshnikov-Kraichnan spectrum in the case of ideal MHD to a Burgers spectrum in the case of a 
shock-dominated hydrodynamic system. The magnetic power spectrum shows a similar behavior. We 
use a ID radiative transfer code to post-process our simulation results; the simulated emission from 
the CS J — 2—1 and H^'^CO+ J = 1 — lines shows that the effects of AD are observable in principle. 
Linewidths of ions are observed to be less than those of neutrals, and we confirm previous suggestions 
that this is due to AD. We show that AD is unlikely to affect the Chandrasekhar- Fermi method for 
inferring field strengths unless the AD is stronger than generally observed. Finally, we present the 
first fully 3D study of the enhancement of AD by turbulence, finding that AD is accelerated by factor 
2-4.5 for non self-gravitating systems with the level of turbulence we consider. 

Subject headings: Magnetic fields — MHD — ISM: magnetic fields — ISM: kinematics and dynamics — 
stars:formation 



1. INTRODUCTION 

The magnetic flux problem is one of the classic prob- 
lems in star formation: How is it that stars form with 
magnetic fluxes many orders of magnitude less than the 
gas out of which they originated? Ambipolar diffusion 
(AD), the process by which the neutral component of a 
magnetized plasma can move relative to the ion compo- 
nent, was identified as a crucial mechanism for resolv- 



ing this problem m ore than half a century ago (Mes 



tel & Spitzer 1956). Subsequent study of the quasi- 
static AD-driven~collapse of an isolated cloud core (e.g 



Spitzer|[T968 l |Mouschovias| |T976l [T9771 [T9791 [Nakano fc 
Nakamura|| 1978; Shu 1983; Fiedler fc Mouscliovias||1992| 



1993p showed that AD could indeed reduce the flux by 
a large factor. However, both observations (e.g. |Zuck-' 



erman & Evans '1974; 'Zuckerman & Palmer"19741 and 
theory ("eTg. Arons fc Max 1975; Mac Low & Klesseli| 
2004; Ballcstcros-P aredes et al. 2007; McKce & Ostriker 
2007) have revealed that supersonic turbulence in molec- 
ular clouds plays a central role on the whole star for- 
mation proces s. Since molecular clouds are magnetized 



(e.g. IC rutchcr 11999) [Heiles k Troland||2005[ [Troland fc 
[Crutchcr 2008), we are confronted witn the probfem of 
understanding AD in a turbulent medium. 
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While important insights have been obtained through 
analytic studies of AD in a turbulent medium in two di- 
mensionsJZwcibcl 2002; Fatuzzo k Adam s 2002; Kim fc| 
Diamond]2002; Hcitsch et al.^200l |, the only feasibfe ap- 
proach for addressing this probfem in three dimensions 
is numerical simulation. Unfortunately, even numerical 
simulation faces difficulties in the study of turbulent, 
magnetized molecular clouds with AD: If the neutrals 
and ions are treated as a single fiuid, then AD becomes a 
nonlinear diffusion problem with a time step that scales 
as the square of the grid spacing, which becomes increas- 
ingly challenging at high resolution. On the other hand, 
if the neutrals and ions are treated as separate fluids, 
then the time step scales i nversely as the fastest signal 
velocity . Observations (e.g. Caselli et al.|1998 Bergin et 



al. 1999( ) show that the ionization fractioii in star-forming 



regions of molecular clouds is typically < 10^^. For such 
low ionization fractions, the ion Alfven speed is at least 
three orders of magnitude larger than the neutral Alfven 
velocity, which leads to a time step that is smaller than in 
the hydrodynamic case by the same factor. A two-fiuid 
simulation of AD in a star-forming region would require 
a time step > 10'^ times smaller than for an ideal MHD 
simulation; on a 512^ grid, a simulation lasting several 
dynamical times would require ^ \{f CPU hours, even 
on the most advanced parallel computers using the most 
efficient massively parallelized MHD codes. 



2 



Li, McKee, & Klein 



Li, McKee, fc: Kleiri| ( 2006 1 developed the Heavy-Ion 
Approximation tor MHlJ witii AD and tested it with the 
ZEUS-MPA D code (this approximatio n was used inde- 
pendently by [Oishr^Mac^Low] (|200^ The Heavy-Ion 
Approximation reduces the ion Altven velocity by in- 
creasing the mass of the ions, while at the same time 
reducing the coupling to the neutrals so as to maintain 
the same rate of momentum transfer between the two 
species. This approximation is valid onl y when the ion 
inert i a is n egligible. As pointed out by |Tilley fc Bal-| 



sara 



( 2008 1 , it is therefore impossible to represent high 

trequency MHD waves that propagate in the ion fluid 
with this approximation. However, this is a positive fea- 
ture of our method, not a flaw: it enables simulations 
of very weakly ionized plasmas, such as those in molec- 
ular clouds, to proceed without having to follow waves 
with phase velocities that can reach thousands of km 
in a fluid with v elocities o f only km s~^. Contrary 
to the statement by Li et al. (20101 that "ion-neutral 
decoupling does not appear m the two-fluid treatment 
with heavy ion approximation," our series of papers on 
turbulent AD clearly demonstrates that ion-neutral de- 
coupling clearly appears in two-fluid treatment when the 

Heavy Ion Approximation is used correctly^ 

Through a series of 256^ simulations, Li et al. (20081 



(hereafter Paper I) showed that this approximation can 
successfully accelerate sub-Alfvenic turbulence simula- 
tions with AD by a factor of 100. They studied the 
statistics of non-ideal MHD turbulence with AD using 
two fluids, neutrals and ions. The purpose of Paper I 
was twofold: (1) to study the convergence behavior of 
the Heavy Ion Approximation in order to determine the 
optimal increase in the ion mass; and (2) to investigate 
the effect of AD on the statistical properties of a super- 
sonically turbulent system. The results were presented 
as a function of the AD Reynolds number, which is the 
ratio of the characteristic AD timescale to the flow time, 
or, equivalently, the ratio of the size of the system to the 
characteristic AD lengthscale (Zweibel & Brandenburg 



[1997, ), 



i?AD(4) 



(1) 



where 7ad is the ion-neutral coupling parameter, p^, „ is 
the mean density of the ions and neutrals, respectively, 
^0 is the size of the system (for the simulation, it is the 
size of the turbulent box), v is the rms velocity, and 
-Brms is the rms magnetic field. The importance of AD 
increases as i?AD decreases. In Paper I we found that 
the power spectra of the velocity and the magnetic field 
change from those for ideal MHD to those for hydrody- 
namics as i?AD(^o) goes from very large values to values 
less than unity. Most of the simulations in Paper I used 
a unigrid 256'^ in size, with only one 512"^ model to show 
that the results were converged. 

Since stars form from dense clumps, the properties of 
clumps in a turbulent medium are important. To study 
dense clumps that arise in supersonic turbulence in the 
presence of AD, we carried out fo ur more 512'^ models 
with different values of i?AD(^o) in McKee et al. (20101 



(hereafter Paper II). The AD models presented in Pa- 
per II span four orders of magnitude in i?AD(-^o)- It 
should be noted that neither Paper I nor Paper II in- 
cluded gravity. In paper II, we used the CLUMPFIND 



algorithm (Wilhams, De Geus, & Blitz 19941 to deter- 
mine the clumps in our simulations. We detiiied clumps 
as connected regions with a density larger than the mean 
density of the turbulent box (see paper II on details of 
how we used CLUMPFIND to extract clumps from the 
simulations). We found that AD affected the clump mass 
function and mass-to-flux ratio, even in the absence of 
gravity. We introduced the categorization of AD into 
five regimes by comparing the gas flow time, tf — i^/v, 
with the ion-neutral and neutral-ion collision times, 
and tni, respectively. The typical value of i?An( ^n) for 15 
cloud clumps with measured magnetic fields (Crutcher 
1999 1 was found to be ~ 20, in excellent agreement 
with the theoretically expected value for self-gravitating 
clumps. Paper II also presented a detailed discussion 
of scaling relations in simulations with AD and showed 
that simulations that include self-gravity cannot vary the 
self-gravity and AD independently unless the ionization 
is treated as a free parameter. The results of Papers 
I and II show that AD can have important effects on 
the clumps even before the core collapse phase that is 
the focus of most quasi-static AD-driven star formation 
models. 

In work that i s complementary to that in t he present 
series of papers, Kudoh & Basu (2008 2011) have sim- 
ulated MHD turbulence, including AL), m a thin, self- 
gravitating slab of gas. The turbulence was initialized in 
the plane of the slab and allowed to decay. The focus 
of their work was on the time required for bound cores 
to form in gas in which the magnetic field initially domi- 
nated gravity. By contrast, our work considers the effects 
of AD on MHD turbulence in the absence of gravity. The 
turbulence is driven continually in all three dimensions 
so as to maintain a constant Mach number. 

In this paper, we present our results on the observa- 
tional implications of AD in turbulent molecular clouds. 
In §2, we summarize our models and the assumptions 
that we have adopted. In §3, we use our higher res- 
olution models to revisit the power spectra of velocity 
and magnetic field found in Paper I. In §4, we study the 
linewidth-size relation of the two fiuids using different 
statistical methods and compare with other work using 
ideal MHD simulations. §5 discusses the ratio of the line 
widths of the ions an d neutrals, which was suggested as 
an indicator of AD by Houde et al. ( 2000a|b ). One of the 
challenges facing observers studying molecular clouds is 
determining the magnetic field strength inside molecu- 
lar clouds. A commonly used method for doing this is 



the Chandrasekhar - Fermi method (CF method) [Chan- 
drasekhar & Fermi (1953), and we discuss the effects of 
AD on this method in §6. Turbulent enhancement of dif- 
fusion processes is commonly observed in experiments, 
observations, and numerical simulations, and in §7 we 
discuss the turbulent enhancement of AD in molecular 
clouds. We summarize our findings in §8. 

2. MODELS AND ASSUMPTIONS 

The results presented and discussed in this paper are 
from the same models presented in Paper II. There are a 
total five AD models with i?AD incremented by a factor 
of 10 from 0.12 to 1200 using a 512'^ grid. Simulations 
of an ideal MHD model, corresponding to i?AD = co, 
and a pure hydro model, corresponding to i?AD = are 
also included for comparison (see Table 1). All models 
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are dri ven by a fixed pa ttern turbulence field, generated 



by the Mac Low 
Mach number A' 



( 1999 1 recipe, to maintain a constant 
37"^ The turbulence is driven at the 
largest scale between wavenumbers 1 < fc < 2. (The 
dimensionless wavenumber k is related to the physical 
wave number by fc = fcpijys^o/(27r)). The inertial range 
of the turbulence extends from fc 3 to fc ^ 20, where 
the upper limit is set by numerical diffusion. Boundaries 
are periodic and the turbulent region is initially threaded 
by a uniform magnetic field with a plasma /? = 0.1. The 
AD models are based on the assumption of ion conserva- 
tion. Discussion of this assumption, including compari- 
son with that of ionization equilibrium, can be found in 
the Appendix of Paper I and in Paper II. For our AD 
models, there is no clear difference between using ion 
conservation and using ionization equilibrium. Gravity 
is not included in any of the simulations. 

Molecular gas is observed to be approximately isother- 
mal, and we assume that here. In fact, temperature fluc- 
tuations are expected due to the intermittency of turbu- 
lence; in supersonic turbulence, gas is heated in shocks, 
for example. Since the cooling li nes in molecular gas are 
often optically thick (Goldsmith & Langer 1978), accu- 
rate calculation of the temperature is a complex problem 
in radiative trans fer that is beyond the scope of this pa- 
per. We note that Pan & Padoan ( 2009 ) have shown that 
cosmic ray heating dominates in most of the volume of 
a supersonically turbulent medium such as we consider, 
and the isothermal assumption holds there. AD heating 
dominates in a small fraction of the gas, however, and 
that will be the subject of a future paper. 

In Paper II, we gave a detailed discussion on how to 
convert the results of simulations to physical units for 
easy comparison to observations, and we do not repeat 
that here. As a reminder, for the AD model with Rab = 
12, the size of the simulation region is 0.41 pc, the total 
mass is 92Mo, the turbulent flow time across the region 
is 0.71 Myr, and the initial mean magnetic field is 90 /iG. 

3. HIGH RESOLUTION SIMULATION POWER SPECTRA 
AND THE CHANGE OF CLOUD MORPHOLOGY 

In Paper I, we discussed the turbulent power spectral 
indexes as functions of Rad based on 256^ simulations, 
although we presented one 512"^ model as part of a con- 
vergence study. We considered a range of Rad{^o) from 
0.12 to 1200. With the doubled resolution of the simu- 
lations in this paper, the inertial range is well developed 
up to at least fc = 20 (Figure [T]) . The resulting spectral 
indexes for the five models using fitting from fc = 3 
to 20 are given in Table 2. All the fitting results and the 
plots shown in this paper are time-averaged results from 
14 data dumps over two crossing times unless specified 
otherwise. The spectral index, n^, for a variable x is 
defined by 

P{k)dk (X k-'^'^dk, (2) 

where P{k)dk is the power in the wavelength interval 
dk. The uncertainties in the power-law fits to the iner- 
tial range slopes are improved. Comparison of the re- 
sults listed in Table 3 of Paper I with the new, higher- 
resolution models shows that most of the results from the 
256^ models are within the uncertainties. This confirms 
the conclusions in Paper I on how the power-law indexes 
of velocity and magnetic field change from ideal MHD 



at very large i?AD(^o) to pure hydrodynamics at small 
^ad(^o)j where ambipolar diffusion dominates. We plot 
the indexes of the ion and neutral velocity power spec- 
tra and of the magnetic field power spectra in Figure 
[2] for better visualization of the transition of the turbu- 
lent systems. At large Rad, both the ion and neutral 
velocity spectra are close to the Iroshnikov-Kraichnan 
(Iroshnikov 1963 Kraichnaii||1965 ) spectrum, which has 
Uy = 3/2. As Rad decreases, the neutral velocity spec- 
trum evolves towards the shock-dominated Burgers spec- 
tnmi (Burgers 1974), which has = 2 and is expected 
for purely hydrodynamic supersonic turbulence. Note 
that for supersonic, super- Alfvenic turbulence, the veloc- 
ity spectral index is a lso similar to a Bur gers spectrum, 
with Uyn = 1-9 ~ 2.0 ( [Padoan et al]|2007[ ). The spectral 
index for the ion velocity also becoines steeper but does 
not appear to be a Burgers spectrum yet, with n^i ^ 1.85 
over the range i?AD(^o) ~ 0.1 — 10. Apparently, the in- 
frequent collisions with neutrals at such small values of 
Rad are not able to turn the ion velocity spectrum into 
a Burgers spectrum. The transition between ideal MHD 
and the AD-dominated, quasi-hydrodynamic behavior of 
the turbulence occurs over the range i?AD(^o) '-^ 1 — 100, 
i.e., at relativel y large AD Rey nolds numbers. In the 

on the propaga- 
imensions, the ef- 



recent study by van Loo et al. 
tion of nonlinear MHD waves in two 
fects of ambipolar diffusion also set in at significantly 
larger scales than one would naively expect at around 

RAD{i) - 1. 

In Paper II, we categorized the effect of AD into five 
regimes: 

I. Ideal MHD {Rad — > oo): The ions and neutrals 
are perfectly coupled. 

II. Standard AD (i/ > t„i '> Un, corresponding to 
Rad > A^A^): The neutrals and ions are coupled 
together over a flow time. 

HI. Strong AD {t^i > tj > tin, corresponding to 

Ma^ > Rad > ^Ai^)' The neutrals are no longer 
coupled to the ions in a flow time, but the ions 
remain coupled to the neutrals. 

IV. Weakly coupled {tin > tf, corresponding to 
A4Ai > Rad)- The ions and neutrals are only 
weakly coupled and act almost independently. 

V. Hydrodynamics {Rad ^ 0): The neutrals are not 
affected by the trace ions and act purely hydrody- 
namically. 

Here tf is the flow time, or the dynamical time scale of 
the turbulence, tni is the neutral-ion collision time scale, 
and tin is the ion-neutral collision time scale. From Table 
2 and Figure [2] the transition between ideal MHD and 
the AD-dominated, quasi-hydrodynamic behavior of the 
turbulence occurs in the Standard AD regime. Near the 
beginning of the Strong AD regime (i.e., for Rad some- 
what less than A4a ), the neutral power spectrum is very 
much the same as a Burgers spectrum. F ffteen molecular 
clumps with measured magnetic fields (Crutcher 19991 
have Rad from 3 to ~ 70, with a logarithmic mean of 
the AD Reynolds number ^ 20 (see Paper II). All these 
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clumps are inside the Standard AD regime, and measure- 
ments of their velocity power spectra should show a tran- 
sition from a spectrum close to the Iroshnikov-Kraichnan 
spectrum to one close to a Burgers spectrum as i?AD(^o) 
goes from the highest observed value to the lowest. The 
power index of the magnetic field goes through a similar 
steepening as AD becomes stronger in regime II, mimick- 
ing closely the transition in the velocity power indexes. 
Previous studies of the density power spectrum f or su 



personic hy drodyn amic and M HD turbulence (e.g. Kim 
fc Ryu|2005 . Kowal et al.|2007| ) have found that the^eu 
sity spectral mdex, Up, is dependent on the sonic Mach 
number but (for MHD turbulence) not very sensitive to 
the Alfven Mach num ber. Using 512^ sup ersonic hydro- 
dynamic simulations, IKim & Ryu (20051 reported that 



Up varies from 0.52 to 1.73 as the Mach number varies 
from 12 to 1.2. In the limit of high Mach number, shocks 
compress the gas into infinitely thin sheets, leading to 
Pp{k) cx pressure forces become important as the 
Mach number decreases, leading to rio ~ 2 at ~ 1 

MHD super- 



(Saichev & Woiczynski 


19961 


From 256^ 


sonic turbulence simulations. 


Kowal et al. 



that Up is shallower than in hydrodynamic supersonic 
turbulence: For ^ 0.7, they found Up = 0.5 ± 0.1 at 
a sonic Mach number M = 7.0 ± 0.3 and Up — 1.3 ± 0.2 
at = 2.2 ±0.03. 

The density spectral indexes for all our models are tab- 
ulated in Table 2. The Alfven Mach number of all the 
turbulence models reported in this paper is Ma = 0.67 
and the sonic Mach number is = 3. First we check 
to see if our models at the extremes of ideal MHD and 
pure hydrodynamic turbulence are consistent with other 
simula tion results. From the results of IKowal et al.l 
( |2007[ ) and |Kim fc Ryill ( |2005l) , we expect Up ~ 1.17 
and 1.2, respecHvelyTlor ideal MHD (model m3i) using 
a simple linear interpolation; consistent with this, we 
find Up = 1.08 ± 0.07. Our AD model m3c2r-l, with 
^ad(^o) = 0.12, should be very close to the hydrody- 



namic limit; we find n„ 



1.23 ± 0.05, consistent with 



the hydrodynamic results of |Kim fc Ryu (|2005|). The 
fact that Uo and nn n for models m3i and in3c2r-l are 



numerically about the same is a coincidence for our par- 
ticular choice of sonic and Alfven Mach numbers. How- 
ever, the amplitudes of the density power spectra for the 
m3i and m3c2r-l models are not the same because AD 
allows larger density contrasts (larger dispersion in the 
density PDF, see Paper I) to occur in the neutral gas 
for model m3c2r-l (i?AD(^o) = 0.12) than in the gas 
for model m3i (ideal MHD). This is most likely because 
compression in the ideal MHD case is primarily along the 
field lines, whereas for low i?AD(^o) compression can 
occur in all three dimensions. 

From Table 2, we see that AD has an interesting effect 
on the density power spectrum as i?AD(^o) varies, even 
though the values of Up for the ideal MHD and pure hy- 
dro regimes are roughly the same for our choice of tur- 
bulence initial conditions. The index Hp steepens from 
1.08 to 1.77 as i?AD(^o) drops from c» to 12 and then 
changes back to 1.23 at i?AD(^o) = 0.12. We attribute 
the steepening between ideal MHD and i?AD(^o) = 12 
to the fact that at high values of i?AD(^o), AD is most 
effective on small scales, and this leads to damping of 
the density fiuctuations on such scales. As -Rad(^o) de- 



creases further, AD becomes effective on all scales and 
the results approach the hydrodynamic limit. These re- 
sults are shown graphically in Figure|3] where we plot the 
density power spectrum of model m3i and the neutral 
density power spectra of models m3c2r2, m3c2rl, and 
m3c2r-l. The neutral density power spectrum of model 
m3c2r3 (i?AD(^o) = 1200) is similar to the density spec- 
trum of m3i, and that of model m3c2r0 (i?AD(^o) = 1-2) 
is similar to that of m3c2r-l, so they are not plotted. The 
range i?AD(^o) = c» to 12 includes the transition from 
AD regime I to AD regime II. 

AD also has a significant effect on the morphology of 
the gas. In Figure |4j we plot density slices from mod- 
els m3c2r3 (i?AD(4) = 1200), m3c2rl (i?AD(4) = 12), 
and m3c2r-l (i?AD(^o) = 0.12) in the x — z plane at 
the same value of y; the mean field is in the z-direction. 
The density slices are at taken at the end of the simula- 
tions and displayed on a linear scale in order to better 
show the changes in the spatial distribution of the high 
density gas. Figure [4^ shows clear filamentary struc- 
tures aligned roughly along the mean i?-field direction 
in model m3c2r3. The curved filamentary structures 
trace the perturbed _B-field lines. In this model, gas can 
freely move along the i?-field but cannot easily cross it. 
When AD becomes important, the morphology changes 
because of diffusion across the i?-field (see Figure l4jD, 
^ad(^o) = 12). The density spectral index steepens as 
discussed above and the maximum density also increases. 
In going from -Rad(^o) = 12 to 0.12, the system moves 
from AD regime II into regime HI, and the morphology 
changes again (see Figure ]4t). Fragmentation becomes 
important because the magnetic field can no longer sup- 
press motions normal to the mean field direction, and 
as a result, high-density structures form at much smaller 
scales and the density power spectrum flattens. The in- 
crease in fragmentation from -Rad(^o) — 12 to 0.12 is 
also reflected in the clump mass function, size, and ge- 
ometry, as discussed in Paper II. The changes in the den- 
sity power spectrum as AD becomes more important are 
thus mirrored in the morphology of the gas. 

Figure [5] shows what observers would see if they ob- 
served the gas in our simulations with the Atacama Large 
Millimeter Array (ALMA). We used the CAS A simdata2 
software to create the simulated observations. CASA 
(Common Astronomy Software Application) is a suite of 
software tools for calibration and analysis of radio as- 
tronomical data (Jaeger 2008) . Simdata2 turns a model 
of the sky into tne visibilities that would be measured 
with ALMA and then produces a synthetic image from 
the model visibilities. The panels on the right column 
are simulated images from models m3c2r3, m3c2rl, and 
m3c2r-l assuming a 2-hour observation of CS J = 2 — 1 
emission at a frequency of 98 GHz. The panels on the 
left are simulated images of the corresponding ion com- 
ponent from the simulation, H^'^CO+ J = 1 — emission 
at frequency of 86.6 GHz, also for 2-hour observation. 
The viewing direction is normal to the mean field direc- 
tion, which is along the vertical axis . We use the ID 



radiative transfer code SimLine ( Ossenkopf et al.| 2001| 
to compute the line emission from a region of a molec- 
ular cloud with size 0.41 pc and density n(H) = 10^ 
cm~^. Since Rab is proportional to the ionization frac- 
tion, the the two limiting cases m3c2r3 and m3c2r-l, with 
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i?AD = 1200 and 0.12 respectively, are models with the 
highest and lowest ionization fractions. The abundances 
of some molecules, such as CS and HCO+, are found 
to be insensitive to high ionization environments, both 
from observations and modeling, even over a two-order- 
of- magnitude chang e in the ionization rate p^arquhar~et 
al.| |1994| |Lepp fc Da lgarno 19961 jMeijerink et al.||2Ulip . 
Modeling of the abundances of CS and HCU+ shows that 
there could even be an anti-correlation between abun- 
dance and ionization rate at high column density. For 
simplicity, we fix t he CS abundance to be 10^^ for all 
models (e.g. Furuya et al.||2011[ ). We assume the mean 
abundance HCO"*" relative to Hq is 2 x 10~^ for the stan- 



dard ionization model, m3c2rl (e.g. Padoan et al.|[2004 
Furuya et al. 2011); this model has ii^AD = 12, close 
to the median value in Crutc her's sample of molecular 
cloud clumps ( Crutcher||1999 and see paper II). For the 
low-ionization model, mc3cr-i, we assume that the mean 
abundance of HCO+ tracks the total ionization, which is 
100 times less. For the high-ionization model, m3c2r3, 
we adopt a mean abundance for HCO+ of 5 x 10~^. We 
assume that th e abundance of H^- ^ CO+) is 2% of tha t 
for HCO+ (e.g. |Purcell et al.||2006| [Furuya et al.||2011[ ). 
Since the ion density varies at ditterent location in the 
turbulent box, the local abundance of H^^CO"'" relative 
to H2 will depend on the local ion density from our sim- 
ulations. The rele vant molecular data are taken from 
the JPL catalogue dPickett et al.| |1998|) and the CDMS 
database (Miill er et al.|200T ). We compute the line emis- 
sion intensities along the line of sight as functions of fre- 
quency from these models using SimLine. Both lines are 
mostly optically thin and have maximum optical depths 
of order unity. The results are used as input for CASA 
simdata2 to generate the images that would be observed 
with ALMA. The size of the turbulent box is assumed to 
be 1' in the maps. When AD is negligible, as in the case 
of model m3c2r3, the simulated images of the ion and 
neutral gases are basically identical, as expected. When 
AD becomes important, the simulated images show that 
the ion and neutral emission differ from each other. The 
simulated images in Figure [5] demonstrate that decou- 
pling of neutrals from ions is observable in principle at 
values of Rab typical of clumps with measurable mag- 
netic fields. Observation of a difference between ion and 
neutral images is not in itself sufficient to conclude that 
AD is operating, however, since the difference could be 
due to chemical effects as well.; if the images are the 
same, however, AD is unlikely to be important. 

4. LINEWIDTH-SIZE RELATION 

Larson (1981) was the first to point out the correlation 
of size and velocity dispersion in molecular clouds, which 
he attributed to turbulence: 



(3) 



He found q = 0.38, but subsequent observations using 



various techniques (e.g., iSolom on et al, 
^aliylfigMl [Lazarian fc ' PogosyaH piMf 
2006f ha^ found q - 0.5. Since 



1987 
Ussen 



Miesch fc 
opt et a 



Mir= I Pv{k)dk, 
lk=io/e 



(4) 



the spectral index of the velocity power spectrum, n^,, is 
related to the linewidth-size parameter q by 



2g + l. 



(5) 



Passot, Pouquet, fc Woodward (|1988) pointed out that 
with q = U.5, the linewidth-size relation corresponds to 
the Burgers spectrum, Pv{k) = k^^. 

Note that the linewidth-size index q in ccjuation ([3| 
is not identical to the one inferred by observers from a 
Princ ipal Components Ana lysis (PCA: |Heyer fc SchloeTb] 



1997 



relation as 



Brunt fc Heyer|2002 ) , who write the linewidth-size 



^v, obs 



cx R". 



(6) 



As discussed by Brunt fc Heyer (2002), the information 
retrieved from a principal component analysis of the ob- 
served spectra involves a line-of-sight projection of the 
emissivity (modulated by the opacity for optically thick 
lines), which depends on the velocity, density, tempera- 
ture, and abundance, all of which vary al ong the line of 
sigh t. (By contrast , the d ata analyzed by Larson ( 1981 1 
and Solomon et al. (1987), for example, focused primar- 
ily on entire clouds, where the sizes are determined by a 



transit ion from molecular to atomic C.) Brunt & Heyer 
(2002) use numerical models to find a conversion be 



tween the observed linewidth-size index a and the veloc- 
ity power spectral index q. For example, for a log-normal 
density distribution, they find that an observed value of 
a ~ 0.6 corresponds to q ~ 0.5. 

Here we compare the linewidth-size relation deter- 
mined using the power spectrum in Fourier space with 
that in physical space. We shall not discuss the alterna- 
tive method of determining the velocity po wer spectral 



index from spectral maps of emission lines ( pLlazarian fc 
|Pogosyan 2000) , which has b een verified by high resolu- 



tion simulations (Padoan et al. 2006) 



In the bo x-decomposition method (e.g. Lemaster fc 
Stone 2009[ ), we divide the turbulent box into succes- 
sively smaller boxes of equal size, £ = £o/2"^, until the 
boxes are of size 2x2x2 grid cells. We use volume 
weighting to compute the velocity dispersions of the ions 
and neutrals in the local frame of each box so that we 
can compute the velocity power spectral index directly. 
Lemaster fc Stone ( 2009 ) also used the volume- weighted 
velo city dispersion in their box deco mposition analysis. 



but Falceta-Gongalves et al. (20101 used the density 



weighted velocity, pv, in their analysis 

In Figure [6j we plot the ion and neutral velocity dis- 
persions of models m3c2r-l (i?AD(^o) = 0.12), m3c2rl 
(i?AD(^o) = 12), and m3i (ideal MHD) as functions of 
length scale based on the box decomposition method. 
Note that in model m3c2r-l, which has strong AD, the 
ion velocity dispersion is smaller than the neutral veloc- 
ity dispersion at all length scales. For model m3c2rl, 
the ion velocity dispersion is just slightly smaller than 
the neutral velocity dispersion (see Section p] for further 
discussion of the ion-neutral line ratios). The ion and 
neutral velocity dispersions for model m3c2r3 are nearly 
the same as those for the ideal MHD model, m3i. 

A well defined power-law region can be identified in the 
linewidth-size relation plotted in Figure [6j. This power- 
law region almost exactly overlaps the inertial range in 
the velocity power spectrum, fc ~ 3 — 20. The time- 
averaged power indexes in the inertial range based on 
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the box decomposition method, qhd, are listed in Table 
3. These indexes are derived from a fit of the three 
data points corresponding to £/£o = 1/4, 1/8, and 1/16. 
This linewidth-size parameter in physical space, gbd, can 
be compared to the one implied by the velocity power 
spectrum in Fourier space, = [uy — l)/2. The values 
of Qps for the ions and neutrals, gi^ips and gn,ips, are also 
listed in Table 3. We can immediately see that there is 
a large discrepancy between qps and gbd- The reason is 
simple: qps is based on only the inertial range, whereas 
each box that goes into the determination of (/bd includes 
all scales smaller than the box size; as shown in Figure[6j 
dissipation reduces the velocity dispersion on scales be- 
low the inertial range. Therefore, to obtain a value of q 
from the power spectrum that is comparable to that from 
the box decomposition method, we need to do an inte- 
gration from the actual velocity power spectrum, as in 
equation Q ; we term the resulting power-law index g^ps 
since it includes the dissipative range. Note that the fit 
for 5dps is based on using the same three values oi £/£o in 
equation (3) as were used in determining gbd- Because of 
the much lower power at large k due to numerical dissipa- 
tion, (7(jps is steeper than ^pg, which is estimated from just 
the inertial range. The results for both ions and neutrals 
are listed in Table 3, which shows that qdps matches q^d 
much better than qps does. From this exercise, we can see 
that it is inappropriate to compare the value of q deter- 
mined using box decomposition with observation because 
the uncontrollable numerical dissipation region at large 
k (which some codes confound with a large bottleneck 
effect) significantly affects the linewidth-size index. 

5. ION-NEUTRAL LINEWIDTH RATIO 
Observations of lines of molecular ions (e.g. HCO+, 
HCS+) and of neutral molecules (e.g. HCN, H^^CN) 
show that in many cases the neutral linewidth s exceed 
the ion linewidths ( |Houde et a"L]|2000a|b[ |2002[ ). This is 
a striking result, since HCU"'" has a lower critical den- 
sity than HCN and therefore should be more spatially 
extended and corres pondingly have a larger linewidth 
( Houde et al.||2000a |. These authors suggested that AD 
is responsible tor this effect and presented a microphys- 
ical model to account for it. They began by assuming 
that they were in the frame in which the electric field 
vanishes. For molecular clouds, in which the ion gyrofre- 
quency is large compared to the collision frequency, this 
corresponds to the ion frame. They showed that if the 
neutrals drift relative to the ions, the ions will acquire 
motions perpendicular to the field that correspond to 
the temperature associated with the neutral drift veloc- 
ity, rriiaf ~ ninv'^j.^f^^. In fact, t his ion heatin g is well 



known in the theory of C-shocks (Draine 1980). Houde 
et al. (2000a) then concluded that because the velocity 
dispersion of the hot ions is less than the drift velocity 
of the neutrals (since the ion mass is much larger than 
the mean neutral mass), the ion linewidths should be 
narrower than those of the neutral molecules. This ar- 
gument completely neglects the turbulent motions of the 
ions (which are associated with hydromagnetic waves in 
the directions normal to the magnetic field ) , how ever, 
and so is invalid. Subsequently, Li fc Houde| ( 2008 ) sug- 
gested a different microphysical model for the ettect of 
AD on the linewidths: they pointed out that the ambipo- 
lar drift velocities needed to account for the observations 



would occur only on small scales, such that Rab ^ 1- 
They assumed that the turbulent velocities of the ions 
would damp out on scales smaller than this, whereas 
those of the neutrals would persist down to the neutral 
viscous scale; as a result, the net velocity disp ersion of 
the ion s would be less than that of the neutrals. ILi et al.l 
(2010) showed that the ion-neutral linewidth ratio was 
smaller in regions of strong magnetic fiel d, which is con 



sistent with AD since i?AD(^o) cx; B~ . Hezareh et al. 
( 2010| assumed that the damping of the ion turbulence 
would occur at a well-defined length scale corresponding 
to i?AD = 1 and inferred the value of t hat length scale in 



DR21 ( 0H). On the theoretical side, Falceta-Gongalves 



et al. (2010) presented arguments supporting the con 



jecture that AD is responsible for the small observed 
values of the ion-neutral linewidth ratio, although their 
simulations, being ideal MHD, could not directly test 
the conjecture. In simulations that are closest to the 



ones presented here, Tilley & Balsara (2010) carried out 
simulations with a two fiuid code and showed that the 
ion-neutral linewidth ratio systematically decreased as 
^ad(^o) decreased, over a similar range as we consider. 
They did not use the Heavy-Ion Approximation, and as 
a result their resolution was limited to 192'^. 
Our results are ful ly consistent with the conjecture by 



Houde et al. (2000a) that AD is responsible for the ob- 
servation that ion linewidths are narrower than neutral 
ones. However, our results do not show the abrupt de- 
cline in the ion velocity power spectrum that they pos- 
tulated, so th ey do not sup port the more detailed model 
put forth by Li & Houde |2008). It must be borne in 
mind that since our simulations are isothermal, we do 
not include the drift heating of the ions, and as a result 
our ion-neutral linewidth ratios are lower limits. 

In Figure [6j we have already seen that the mean ve- 
locity dispersion of the ion component, cr^ i, is always 
less than or equal to the mean velocity dispersion of 
the neutral component, ay^n, at all length scales. The 
stronger the AD effect (smaller Rab), the smaller cr„_i 
is compared to (T„^„. We can plot the ion and neu- 
tral velocity dispersions at different length scales using 
box decompo s ition in the same way as in Figure 3 of 
Houde et al. (2002) for comparison. In Figure [7 a,, we 
plot the volume-weighted mean cr^^i versus the volume- 
weighted mean a^^n computed using box decomposition 
for models m3c2r3, m3c2rl, and m3c2r-l, correspond- 
ing to i?AD(4) = 1200, 12, and 0.12, respectively. The 
error bars show the errors of the time-averaged mean. 
The velocity dispersions from the different length scales 
are plotted: the highest velocity dispersions are from the 
largest boxes {£o /4) , and the lowest are from the smallest 
(£o/128). The highest three points are inside the inertial 
range. We can see that the data points from length scales 
smaller than the inertial range also follow the same trend 
as those inside the inertial range. This result is c onsistent 
with t he observed linewidth ratios reported in |Li et al.] 
( 2010 ). The solid hue corresponds to Oy^ijoy^n = 1- Fig- 
ure^, shows that AD reduces the ion-neutral linewidth 
ratio below unity; for i?AD(^o) — 0.12, which is close to 
the hydrodynamic limit, we find Oy^ijoy^n = 0.62 ± 0.01. 
Figure (Tja shows that this result remains valid when the 
density-weighted velocity dispersions are used instead; in 
this case, Oy^ijoy^n = 0.60 ± 0.01 for i?AD(^o) = 0.12. 
However, observers are not able to obtain 3D veloc- 
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ity dispersions as in the block decomposition method. 



Therefore, we also plot the ratio cr^_; 



of the 



density-weighted velocities normal to the mean field pro- 
jected through the whole box in Figure Up (solid sym- 
bols) . The ratios are the mean values along the two 
cardinal directions normal to the mean field. The ra- 
tios from the projected ay^i and a^^n show the same 
behavior as a function of Rab as the ratios from the 
block decomposition method. It should be noted, how- 
ever, that the situation is quite different for observations 
parallel to the mean field; in that case, the ion-neutral 
li newidth ratio is abo ut unity. Comparison of Figure 3 
in |Houde et al. (|2002 ) with Figure [Tja then suggests that 
most ot the molecular clouds they observed correspond 
to i?AD(^o) ^ 12, which is consistent with the logarith- 
mic mean i?AD(^o) ~ 20 for observed mol e cular clumps 
that we found in Paper II. Houde et al. (2002) found 
some clouds with ion/neutral velocity-dispersion ratios 
considerably below the average, and it is possible that 
these indicate regio ns of strong A D. 



As noted above, Li & Houde (2008) suggested that 
the ion-neutral line width ratio is less than unity be- 
cause the ion velocity dispersion is damped below the 
AD lengt h scale whereas the neutral v elocity dispersion 
is not. |Falceta-Gongalves et al. ( 2010) presented re sults 
that support this idea, and I'illey & Balsara (2010) pre- 
sented direct two-fluid simulations that demonstrate its 
validity. Ou r higher reso lution results are consistent with 
the work of Tillcy & Balsara (2010). However, in contrast 
with the suggestion of ^Li fc Houde (2008J , AD is unlikely 
to produce a sharp break in the ion velocity power spec- 
trum. Even if we assume that the turbulence satisfies a 
line-width size relation and that the gas is in ionization 
equilibrium, Rab is still a function of three free param- 
eters (density, velocity, and magnetic field strength), all 
of which vary in a turbulent medium. Therefore, instead 
of finding a sharp change in the slope of the ion veloc- 
ity power spectrum, we expect a gradual change from 
an Iroshnikov-Kraichnan spectrum to one that is close 
to a Burgers spectrum over ^ 2 orders of magnitude of 
^AD (-^o ), as shown in Figure [2] None of the AD simula- 
tions of [Oishi fc Mac Lo^ (|20D6D [Downes fc O'Sullivan 
( |2009[ ), or |Tilley fc Balsara|poTo| , our 256^ (Paper I) or 
512'^ simulations reveal a change in spectral slope in the 
inertial range of the calculation. One would need an ex- 
tremely high resolution simulation with an inertial range 
spanning 2 ~ 3 orders of magnitude of length scale to 
observe the gradual change in the ion spectrum. 



6. THE EFFECT OF AMBIPOLAR DIFFUSION ON THE 
CHANDRASEKHAR-FERMI METHOD 

Polarimetry is a commonly used method to determine 
the magnetic field orientation in the plane of t he sky 



e:E. 



I'hompson et al. 



Hildebrand et al. 2000 Novak et al. 2000 Ward 



Based on a dynami cal method 



proposed by [Chandrasekhar & Fermi (19531, the mag- 
netic field strength in the plane ot the sky can be esti- 
mated provided that the turbulence is isotropic and that 
the perturbed magnetic energy is in equipartition with 
the turbulent kinetic energy. If one can also measure the 
magnetic field along the line-of-sight using Zeeman map- 
ping, one can obtain the true magnetic field strength and 
the 3D orientation of the magnetic field. 



(2001 1 have carried out a systematic in- 



Heitsch et al. 

vestigation of the accuracy of the Chandrasekhar-Fermi 
(OF) method using 3D simulations of ideal MHD turbu- 
lence. They proposed a new recipe for using the CF 
method to estimate the mean magnetic field strength 
that allows for non-equipartition of turbulent magnetic 
and kinetic energy, large-amplitude fiuctuations of the 
magnetic field, and the limitations of observational reso- 
lution. They define the ratio of the field inferred via the 
CF method to the actual field to be 



acF 



CF 



B) 



(7) 



and conclude that, for strong fields at least, acF — 2 
when the angular structure of the field is well resolved. 

Heitsch et al.'s (2001) result for acF allows for de- 
viations from equipartition between the energy in non- 
thermal motions and that in fluctuating magnetic fields, 
which they parametrize by 



5Up_ 



E 



K 



1 



(8) 



where 5Ub is the energy in the fiuctuating magnetic field 
and Ek is the energy in turbulent motions. Their anal- 
ysis of the CF method is based on the assumption that 
^ is known. The value of ^ has a significant effect on 
the strength of the inferred field: small values of ^ lead 
to small fiuctuations in the direction of the field; if it is 
then assumed that ^ = 1, as in standard applications of 
the CF method, the smallness of the fluctuations is in- 
terpreted being due to a strong field, with acp = ^"-^/^ 
(Heitsch et al. 2001). Both observation and theory are 
consistent with approximat e equipartition (i^ ^ 1). Al- 
though Heitsch et al.| (|2001[) suggested that observations 
by Crutcher (1999) do not support equipartition, the me- 
dian Alfven Mach number of the 15 clouds in Crutcher's 
sample with measured magnetic fields is unity, which is 
consistent with ^ ~ 1 since {SB"^) ~ (_B^) for Ma ~ 1- 
(It should be noted, however, that this simple argument 
allows ^ to deviate from unity by a factor ~ 2.) The- 
oretically, it has been shown that equipartition applie s 
to weak Alfvenic turbulence (IZweibel fc McKeel|1995[), 
including t hat modeled in the theory of Coldreich fc Srid- 



har (19971. However, there is good evid ence for deviaT 
tions from equipartition in simulations. Heitsch et al. 
(2001 1 were agnostic as to whether these deviations are 
a physical effect or an artifact of their simulations. One 
factor that appears to affect the value of ^ is the driv- 
ing scale: Simulations in which the drivii ig is on large 



scales and the field is strong (/3 < 0.2) (e.g . Heitsch etal 



2001] lOstriker et al.||2001[ jlemaster fc S tone 200^1ind 
this work) have ^ < a few tenths; for example, we find 



^ ^ 0.2 in our ideal MHD turbulence model 
the other hand 



m3i. On 



(1998) 



the simulations of Stone et al 
in which the driving was on a scale small compared to 
the box scale, found ^ ~ 0.6 for /3 < 0.1. This result can 
be understood as being due to approximate equipartition 
between the magnetic fluctuations and the turbulent mo- 
tions normal to the field lines; for isotropic turbulence, 
this would lead to ^ ~ |. In our simulations, which are 
driven at large scales, we have found that the reduction 
in ^ appears to be due at least in part to motions in which 
the gas rotates around the field lines. While this may oc- 
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cur in simulations with periodic boundary conditions, it 
is unlikely to occur in Nature. 

Here we address a separate issue: How does AD affect 
the CF method? In Paper I, we showed that AD reduces 
the turbulent magnetic field energy below that in ideal 
MHD turbulence (see Figure 7 in that paper). That is, 
in addition to the physical and numerical effects that 
cause ^ to differ from unity in the ideal case, AD reduces 
^ because the ions are no l onger as well couple d to the 
dominant neutrals. Whereas Heitsch et al. (2001 ) studied 
a range of field strengths, all our models have relatively 
strong fields, with a plasma-/? parameter /3 = 0.1. The 
level of turbulence is modest, with Ai = 3 and A4a — 
0.67. 

Our treatme nt of the CF metho d follows the proce- 
dures given in |Heitsch et al] ( |2001[ ), and we use similar 
notation. In terms ot the Stokes parameters U and Q, 
the polarization angle </> is given by 



1 U 
-arctan— , 

2 g' 



(9) 



When viewed along the y-direction (recall that the mean 
field is in the z-direction) , the polarized intensity is 



P = Q + iU <x / p{y) 



cos^ 7 dy, (10) 



where 7 is the angle between the magnetic field and the 
plane of the sky. The CF estimate of the mean field is 



Bc^^^^)^^ = acAB) 
cr(tand) 



(11) 



where (p) is the mean density, 5 = {(f)— {(/))) is the differ- 
ence between the local and the mean polarization angle, 
and a is the standard deviation. Note that all quantities, 
including the line-of-sight velocity dispersion, a{vios), are 
density- weighted. Because the turbulence is relatively 
weak in our simulations (A^a = 0.67), the perturbed 
field is less than the mean field. To maximize the effect of 
AD, we observe the data perpendicular to the mean-field 
direction, computing U and Q from B^ and Bz when ob- 
served in the y-direction (and By and B^ when observed 
in the a:-direction) . We take Bqf to be the mean of the 
estimated field strengths along the x and y directions. 

In the presence of AD, the departure from equiparti- 
tion will be larger than in the ideal case because of the in- 
efficient transfer of turbulent kinetic energy to magnetic 
energy, as wa s demonstrated analytically by |Zweibel"^ 
McKee ( 1995 1 and numerically in Paper I (see Figure 7). 
In order to tocus on the effect of AD without the com- 
plications introduced by deviations from equipartition in 
the ideal MHD case, we normalize the values of ocf from 
our simulations by the value of qcf from the ideal MHD 
model m 3i. (For this model, ac F — 2, similar to the value 
found by [Heitsch et al.| ( |2001[ ).) The results in Figure [s] 
show that AD starts to affect the CF method somewhere 
in the range i?AD(^o) = 1 ~ 10. At i?AD(^o) ^ 1, we 
find that acF ~ 3 times the ideal MHD value. In other 
words, if we assume that our estimate of the CF coef- 
ficient, acF, is correct for the ideal case, then applica- 
tion of the CF method to a system with i?AD(-^o) — 1 
would give a field strength that is too large by a factor 
3. At i?AD(^o) = 0.12, the normalized CF coefficient is 



Zweibel 


( 1990 1 also 


;h wouk 


be overes- 



acF > 17. Thus, a large correction 
will be required if AD is very strong, 
found that the magnetic field streng 
timated in the presence of AD using a highly simplified 
clump and intcrclump gas model. However, the 15 molec- 
ular clouds with measured magnetic fields discussed in 
Paper II have a logarithmic mean i?AD(^o) ~ 20, and 
most of them have i?AD(^o) > 10. According to Fig- 
ure [Sj AD has only a minimal effect on the inferred field 
strength, Bqf, for such clouds, and it is safe to estimate 
the mean B field based on ideal MHD models. For re- 
gions where i?AD(^o) 1, a significant correction for AD 
would be necessary. 

7. TURBULENT ENHANCEMENT OF AMBIPOLAR 
DIFFUSION 

Enhancement of diffusion processes by turbulence is 
observed in ex periment, obse rvations, and numerical sim- 
ulations (e.g. JCowee et al. 



2009 Liu et al. 2009 Na- 



gai et al. 2009 H 'I'he enhancement ot AD by turbu- 



lence h a s been discussed by Zweibel (2002), Fatu zzo fc 
(2002[), |Kim fc Diamondi ( |2002| ), and |HeitscFet 
1). TIach ot these works niade substantial simpli- 



lAdams 

|aL (2004, 

fications: jFatuzzo fc Adams ( 2002| considered field fluc- 
tuations in one dimension, and the remaining authors 
assumed that the motions were confined to two dimen- 



sions; in addition, Kim & Diamond (2002) assumed that 
the gas was incornpressible. Mere, we shall present the 
first fully three-dimensional treatment of the effects of 
turbulence on AD. 

There are two distinct effects of turbulence in the case 
of AD: First, turbulence increases the AD drift velocity, 
Ud, and second, it mixes regions on finer and finer scales, 
so that less diffusion is needed to smooth out variations 
in the mass-to-fiux ratio. The first effect determines the 
rate of AD in a turbulent medium and can be evaluated 
for steadily driven turbulence; that is the objective of 
this section. The second is intrinsically time dependent, 
and depends on the initial conditions; we do not address 
it here, although the effects of turbulent mixing are nat- 
urally included in our simulations. The combination of 
the two effects is turbulent AD, and Zweibel^ (2002) and 



Kim fc Diamond ( 2002 ) have shown that when the neu- 
trals and ions are well coupled (large Rau), the diffusion 
rate for turbulent AD is of order the classical turbulent 
value, v£q, where v is the turbulent velocity on the scale 
4. 

In molecular clouds, the ionization is generally so low 
that the ion inertia is negligible. In that case, the relative 
velocity between the ions and neutrals is determined by 
the balance between the Lorentz force on the ions and 
the drag due to the neutrals. 



jADPiPnVd = 7*- I ( V X B) X B I 

bL, 



(12) 
(13) 



which defines is- Note that this relation is exact in the 
limit in which the ion inertia is negligible, and the defini- 
tion of here differs from the approximate expression in 
previous papers in this series. In terms of the diffusivity 



A; 



B^ 



An JAB Pi Pn 



(14) 
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the drift velocity is 



A 



(15) 



Equation ( 15 ) shows that the effects of turbulence on the 
drift velocity can be divided into two parts, the effect on 
the diffusivity, A, and the effect on the length scale of 
the Lo rentz force, £b- The study by 'Fatuzz o fc Adams| 
( 2002 1 focused on the diffusivity: they pointed out that 
fluctuations in the density and magnetic field in a tur- 
bulent system can increase A and thereby enhance the 
AD rate. Their studies suggest that turbulent fluctua- 
tions can enhance the AD rate by a factor A ^ 1 — 10; 
for mo l ecular clouds, they estim ate d A ^ 2 — 3. |Zweibel| 



( |2002 ), |Kim fc Diamond] ( |2002[ ) and 
considered both the ettect of turbu 



Heitsch et a' 



mm 



ence in reducing Ib 
and the effects of turbulent mixing, albeit in 2D. They 
found that the turbulent AD time is of order the eddy 
turnover time for large R\r)^ which c orresponds to a 
very large enhancement. Zweibel (2002) argued that this 
large enhancement could explain the observed flatness of 
the relation between magnetic fleld strength and number 
density (the B — n relation) in interstellar clouds. 

With the simulations we have performed, we can ad- 
dress the issue of how turbulence affects the rate of am- 
bipolar diffusion through its effects on the drift velocity. 
Our results include both the effects of turbulence on the 
diffusivity, A and on the length scale, is- A direct mea- 
surement of the amount of AD occurring in the turbulent 
box is the value of of the mass flux relative to the mag- 
netic field, PnVd — pvd, averaged over the volume of the 
box. In order to determine the enhancement factor due 
to turbulence, we must have a non-turbulent reference 
value. For large Rad, we define the reference value of 
the drift velocity, Vdo, as the value of the drift velocity 
in a uniform medium in which the field varies on the box 



size, €o- 



VdO 



f>2 

rms 



47r7AD/3iPn4 i?AD(4)' 



(16) 



Note that we have used the rms field strength here so as 
to exclude from the enhancement factor the effects of the 
field amplification that occurs in a high-/? plasma. For 
small values of i?AD(^o): the drift velocity cannot exceed 
the rms velocity, v, so we generalize the expression for 
the reference velocity to 

The enhancement factor for AD due to turbulence is then 
(PnWd) _ H-i?AD(4) / pda 



A = 



PnVdO 



Ra 



PitB 



(18) 



Values of the time-averaged turbulence enhancement 
factor A are listed in Table 4 for all five AD models. 
For i?AD(-^o) > 12, the enhancement is a factor of 2 to 
4.5. For i?AD(^o) < 1-2, the enhancement is negligible 
(A ~ 1): the ions and neutrals are sufficiently decoupled 
that turbulence does not have a significant effect on the 
AD rate (the apparent divergence in A due to the factor 
-Rad(-^o) in the denominator of eq. 18 is more than can- 



celled by the increase in the magnetic scale length, 
gradients in the field cannot be maintained in a weakly 



coupled plasma). Our results for the turbulent enhance- 
ment of AD should apply t o non self-gravi tating regions 
in the ISM, as studied by " " " 
do not apply to AD in sel: 
st ar-forming cores 



Zweibel] (|2002|). Our results 



gravitating regions, such as 



Kudoh & Basu (2008 



201ip have shown that large- 



scaie turbulence accelerates the formation of magneti- 
cally supercritical cores (i.e., cores in which gravity dom- 
inates the magnetic field) by compressing the gas. The 
turbulence was i nitialized in the plane of the slab and 
allowed to decay. Kudoh & Basu (2011) found that the 
time required for botind cores to form scales as the in- 
verse square root of the peak density produced by tur- 
bulent compression, prior to the time that gravity was 
dominant. Since the ambipolar diffusion time and the 
gravitational free-fall time both scale as p^^/^ for the 
typical condition in which the ionization also scales as 
this is just what is expected for a compression in 
the a bsence of turbulence (e.g., Mouschovias & Ciolek] 
1999|). The simulations of [Kudoh fc Basu| (i2008n201 



m 
ekl 



did not address the effect of turbulence on the aD time 
within individual cores, in part because the turbulence 
damped with time. It is possible that the turbulent en- 
hancement factor inside the cores would be reduced by 
the effects of self-gravity since fiux tubes are then not 
as free to move, but confirming this must await a future 
study. 

8. DISCUSSION AND CONCLUSIONS 

In this paper, we have investigated a series of issues on 
the observational implications of molecular clouds from 
the effects of turbulent AD. We summarize our findings 
as follows: 

1. We confirm the results from the 256^ models in Pa- 
per I that the velocity and magnetic power spectra 
are functions of the AD Reynolds number, i?AD- 
The higher resolution 512^ models reported in this 
paper show that the transition of the neutral veloc- 
ity spectrum from the Iroshnikov-Kraichnan spec- 
trum for a strong field, ideal MHD model to the 
Burgers spectrum for a shock-dominated, pure hy- 
drodynamic model occurs primarily in the Stan- 
dard AD regime (Fig. [2]). The ion velocity spec- 
trum diverges from the neutral power spectrum 
when the coupling between ions and neutrals be- 
comes weak, for i?AD(^o) ^ 10. The magnetic field 
power spectrum shows similar behavior. The 15 
molecular clump s with measured magnetic fields 
(Crutcher 19991) are all inside the Standard AD 



regime, and their velocity power spectra should 
show this transition when the high-i?AD(^o) clumps 
are compared with the low-i?AD(^o) ones (see paper 
II). 

2. The density power spectrum steepens as AD 
becomes important and reaches a maximum at 
-Rad(^o) ^ 10; it then becomes shallower as the 
neutrals decouple from the ions at still lower val- 
ues of i?AD(^o)- The steepening occurs because 
AD is most effective on small scales, which leads to 
damping of the density fluctuations on these scales. 
When the neutrals decouple from the ions, neutral 
gas clumps fragment and the spectrum approaches 



10 



Li, McKee, & Klein 



the hydrodynamic limit. The change in the den- 
sity spectral index corresponds to the change in 
the morphology of the gas. We have applied the 
CASA software simdata2 to determine the intensity 
of emission in the CS J = 2-1 and HCO+ J = 1-0 
transitions implied by our simulations. The results 
demonstrate that the decoupling of ions and neu- 
trals is observable in principle. 

3. We have investigated the linewidth-size relation 
from our simulation models using a 3D box decom- 
position method. This method distorts the physical 
linewidth-size relation because the power in a box 
is reduced by numerical damping at small scales. 
The best way to calculate the linewidth-size rela- 
tion in the inertial range is to use equation[5]to infer 
the exponent in the linewidth size relation ixom the 
power-law index of the velocity power spectrum. 

4. Our high-resolutio n two-fluid simulation s confirm 



the suggestion of Houde et al. (2000a) that the 



reduction of the observed lon-neutrai velocity dis 
persion ratio to values below unity is the result of 
AD. However, our results do no t support the micro- 
physical model of Li & Houde (2008) in which the 



damping of the ion motions sets in at a well-defined 
length scale; instead, we find that the effects of 
AD become increasingly important over ~ 2 orders 
of magnitude in length scale. It should be noted 
that our isothermal simulations do not include drift 
heating of the ions, and therefore provide lower lim- 
its on the actual ratio of the ion-neutral velocity 
dispersion ratio when observed normal to the field; 
on the other hand, this field orientation is the opti- 
mum one for observing the effect, and we find that 
the ion-neutral linewidth ratio is about unity when 
observed along the field. Our results also show that 
most of the ob served ion-neutral linewidth ratios 
in the Li et al. ( 2008j sample of molecular clouds 
correspond to /^ad ?^ 10, consistent with the com- 
puted i?AD of the 15 molecular clumps discussed in 
Paper H. 

5. We have examined the effect of AD on the 



Chandrasekhar-Fermi method, which is commonly 
used to estimate the magnetic field strength in the 
plane of the sky from observations of polarization 
angles and velocity dispersions. When AD is signif- 
icant, the magnetic field will be less perturbed and, 
therefore, the corresponding mean field strength 
will be over-estimated. Our analysis shows that AD 
does not significantly affect the CF estimate of the 
field in the range of i?AD of the molecular clumps 
in Crutcher's (1999) sample, There is no need to 
correct values of the magnetic field estimated us- 
ing the CF method for the effects of AD except in 
the very small regions in which -Rad(^o) ^ 1- 

6. With the heavy-ion approximation, we have been 
able to carry out the first fully three-dimensional 
study of the enhancement of AD by turbulence. 
Our simulations show that when the ions and neu- 
trals are reasonably well coupled (i?AD(^o) > 12, 
as is the case for most of the molecular clumps an- 
alyzed in Paper II), AD in non self-gravitating re- 
gions in the ISM is accelerated by a factor ~ 2 — 4.5 
for the conditions we consider. In the weak ion- 
neutral coupling regime, AD is so strong that tur- 
bulence cannot significantly enhance the rate of 
AD. 
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Fig. 1. — Compensated velocity and magnetic field power spectra of 256'^ (dashed line, model m3c2 in Paper I) and 512'^ models with 
-Rad('^o) = 1-2 (solid line, model m3c2r0 in this paper). See Table 3 in Paper I and Table 2 in this pajjcr for the compensated power 
indexes. The inertial range in 512^ resolution extends up to ~ fc = 20. The error bars of the inertial range are plotted in the middle panel 
to show the typical sizes of the error bars in the data fitting. 
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Fig. 2. — Velocity and magnetic field power spectral indexes of five AD models. The neutral velocity power spectral index (squares) 
changes from an Iroshnikov-Kraichnan spectrum to a Burgers spectrum in going from AD regime I to regime III. The ion velocity spectral 
index (circles) also increases as AD becomes important, but it diverges from the neutral velocity index and stays at about n = 1.85 as the 
ions decouple from the neutrals. The magnetic field power spectral index (diamonds) has a similar trend as Rao decreases. 
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Fig. 3. — Time-averaged density power spectrum for ideal MHD (dashed; model m3i) and neutral density power spectra for -Rad(^o) = 120 
(dotted; m3c2r2), -f?AD(^o) = 12 (solid; model m3c2rl), and -Rad(^o) = 0-12 (dot-dashed; model m3c2r-l). The spectral index Up steepens 
to a maximum as -Rad(^o) decreases from infinity (ideal MHD) to -Rad(^o) = 12 and becomes shallower again as -Rad(^o) decreases further. 
See ^for discussion of the relation between the change in the morphology of molecular clouds and the density power spectrum. 
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Fig. 4.— Density slices from models (a) m3c2r3 (Kad(4) = 1200), (b) m3c2rl (i?AD(4) = 12), and (c) m3c2r-l (/?ad(^o) = 0.12), 
showing the morphological changes in the density distribution as the effect of AD becomes stronger. Panel (a) shows sheet and filamentary 
distributions mainly along the mean direction of the magnetic field. Panel (b) shows diffusion of the sheet and filamentary distributions. In 
panel (c), many high-density small structures form as the result of fragmentation. See [jsjfor discussion of the relation between the change 
in the morphology of molecular clouds and the density power spectrum. 
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CS (J=2-l) Hi^CO^ (J=l-0) 




Fig. 5. — Simulated sky images, created from simulation data and post-processed with the radiative transfer code SimLine, that would 
be observed with the ALMA telescope. This figure shows simulated images of the ions (deduced from the H^'^CO"'" J = 1-0 transition — loft 
column) and neutrals (deduced from the CS J=2-l transition — right column) in a molecular cloud with a density of 10** cm""^ as observed 
by ALMA for 2 hours. The top panels are from model m3c2r3 {Rad{^o) = 1200), the middle panels are from model m3c2rl (-Rad(^o) = 12), 
and the bottom panels are from model m3c2r-l (-f?AD(^o) = 0.12). With very weak AD (top row), the images of the ions and the neutrals 
will appear almost the same (^3|. The shapes and sizes of the beams are shown at the left bottom of the images. (Color images available 
in the electronic version.) 
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Fig. 6. — Time-averaged linewidth-size relation for models m3c2r-l (dot-dashed; RadC^o) = 0-12), m3c2rl (dashed; i?AD(^o) = 12), 
and model m3i (solid; ideal MHD). The ion linewidth-size relation is marked by squares and the neutral linewidth-size relation is marked 
by circles. The velocity dispersions are computed using the box decomposition method. The linewidth-size relations obey a power law in 
the inertial range (the region between the two dotted vertical lines) from the power spectra. See fj4] for a discussion of the linewidth-size 
relation for different models. The power-law indexes of the linewidth-size relation for each model are listed in Table 3. 



18 Li, McKee, & Klein 



(a) (b) 




o /c o /c 
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Fig. 7. — Plots of ion and neutral velocity dispersions, (T„ ^ and (T„^„ respectively, for models m3c2r-l (circles), m3c2rl (squares), and 
m3c2r3 (diamonds), (a) Volume-weighted mean (T„ ^ and (T„ „ obtained from the box decomposition method, (b) Density-weighted mean 
cr„ i and cr^.n obtained from the box decomposition method. The solid line corresponds to a^^ija-u^n = 1. The 3 solid symbols are the 
velocity dispersions of the corresponding models projected normal to the mean field. See [jslfor discussion. 
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Fig. 8. — Normalized correction factor acv = Bqy/{B), which is the ratio of the estimated magnetic field strength using the CF method 
to the true magnetic field strength. The correction factor is normalized to the correction factor from the ideal MHD model m3i to show 
the effect of AD on the CF method. The error bars show the uncertainty of the time-averaged mean of the measurements from all the data 
dumps over two crossing times. For -Rad(^o) ^ 10, AD has a negligible effect on the inferred magnetic field strength. AD starts to have 
an effect on the correction factor when 1 < Rad(^o) < 10, and the correction factor is ~ 3 at Rad{^o) = 1- For Rad(^o) < 1, the field 
strength would be significantly overestimated using the CF method. 



TABLE 1 

Model Parameters and Regimes of AD 



Model" 7AD -Rad(^o) Regime of ad 



m3c2r-l 


4 


0.12 


III 


m3c2r0 


40 


1.2 


II ~ III 


m3c2rl 


400 


12 


II 


m3c2r2 


4000 


120 


II 


m3c2r3 


40000 


1200 


I 


m3i 


oo 


oo 


I 



" Models are labeled as "mxcyrn," where x is the thermal Mach number, y = \ logXiol) ^■iid " = log(^AD(^)/l-2). Model "m3i" is an ideal 
MHD. Model m3c2r0 is the same as model m3c2h in paper I. 
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TABLE 2 

Spectral indexes of Velocity and Magnetic Field Power Spectra with Varying Rau 





m3c2r- 


-1 




m3c2r0 




m3c2rl 




m3c2r2 




m3c2r3 




m3i 




RadW)) 




0.12 






1.2 




12.0 




120 




1200 




OD 






1, 


.78±0, 


,06 


1 


.64it0.06 


1 


.SOitO.OS 


1 


.74it0.06 


1 


.54±0.05 


1 


.48±0. 


.06 


'^vi,z(^) 


1, 


.89±0. 


.03 


1, 


.98±0.03 


1 


.88±0.04 


1 


.40±0.06 


1 


.39±0.06 


1 


.38±0. 


.05 


nyi{k) 


1, 


.85±0. 


.03 


1, 


.86±0.03 


1 


.84±0.03 


1 


.58±0.04 


1 


.48±0.05 


1 


.45±0. 


.05 


^vn,r (fe) 


1, 


.97±0. 


.02 


1, 


.92±0.03 


1 


.89±0.04 


1 


.76±0.05 


1 


.55±0.05 








'T'ViijZ (^) 


1, 


.93±0. 


.02 


1, 


.95±0.03 


1 


.88±0.04 


1 


.41±0.06 


1 


.39±0.06 








nvn(fc) 


1 


.96±0. 


.02 


1, 


.94±0.03 


1 


.89±0.03 


1 


.58±0.04 


1 


.48±0.05 








nB,r(fc) 


1, 


.43±0. 


.03 


1, 


.48±0.05 


1 


.45±0.06 


1 


.23±0.05 


1 


.10±0.05 


1 


.14±0. 


.07 


nB,z(fc) 


2, 


.07±0, 


.07 


2 


.14±0.07 


2 


.12±0.05 


1 


.88±0.06 


1 


.64±0.05 


1 


.61±0. 


.08 


nB(fc) 


1, 


.49±0, 


.04 


1 


.56±0.05 


1 


.55±0.06 


1 


.31±0.05 


1 


.17±0.05 


1 


.23±0. 


.07 


npn(*;) 


1, 


.23±0. 


.03 


1, 


.44±0.05 


1 


.77±0.03 


1 


.34±0.03 


1, 


.08±0.06 








npi(fe) 





.83±0. 


.04 





.62±0.04 





.68±0.03 





.78±0.04 


1 


.02±0.05 


1 


.08±0. 


.07 
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TABLE 3 

Line Width-Size Relation Power index q Estimated from Different Methods 



m3c2r-l m3c2r0 m3c2rl m3c2r2 m3c2r3 m3i 



KadW)) 
*,bd" 


0.12 


1.2 


12 


120 


1200 


oo 


0.48 ±0.05 


0.51 ± 0.03 


0.56 ± 0.04 


0.45 ± 0.05 


0.42 ± 0.05 


0.40 ±0.06 


yi,ips 


0.43 ±0.02 


0.43 ±0.02 


0.42 ± 0.02 


0.29 ± 0.02 


0.24 ±0.03 


0.23 ± 0.03 




0.55 ±0.03 


0.57 ± 0.03 


0.64 ± 0.05 


0.52 ± 0.03 


0.46 ± 0.04 


0.43 ± 0.05 



gnbd 0.53 ±0.04 0.55 ± 0.03 0.56 ± 0.03 0.46 ± 0.04 0.42 ± 0.05 
5n'ips 0.48 ±0.01 0.47 ±0.02 0.45 ± 0.02 0.29 ± 0.02 0.24 ± 0.03 
gnjdps 0.60 ±0.04 0.61 ±0.04 0.63 ± 0.04 0.52 ± 0.03 0.46 ± 0.04 



° Estimated from the box-deeomposition method. 
Estimated from the inertial range of the ion velocity power spectrum. 

Estimated from the full range (including dissipation range) of the ion velocity power spectrum. 
The subscript notation also applies to the neutral LWS power index q in the second half of the table. 
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TABLE 4 

Turbulent Enhancement op the AD Rate. 



Model 


Rad(^o) 


A 


m3c2r-l 


0.12 


0.99 ± 0.02 


m3c2r0 


1.2 


1.18 ±0.02 


m3c2rl 


12 


2.23 ±0.09 


m3c2r2 


120 


3.11 ±0.24 


m3c2r3 


1200 


4.53 ± 0.36 



